sep_smooths_year2

Load all required libraries.

library(tidyverse)

## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --

## v ggplot2 3.3.5     v purrr   0.3.4
## v tibble  3.1.6     v dplyr   1.0.8
## v tidyr   1.2.0     v stringr 1.4.0
## v readr   2.1.2     v forcats 0.5.1

## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(plotly)

## 
## Attaching package: 'plotly'

## The following object is masked from 'package:ggplot2':
## 
##     last_plot

## The following object is masked from 'package:stats':
## 
##     filter

## The following object is masked from 'package:graphics':
## 
##     layout

library(broom)

Read in raw data from RDS.

raw_data <- readRDS("./year2.RDS")

Make a few small modifications to names and data for visualizations.

final_data <- raw_data %>% mutate(log_copy_per_L = log10(mean_copy_num_L)) %>%
  rename(Facility = wrf) %>%
  mutate(Facility = recode(Facility, 
                           "NO" = "WRF A",
                           "MI" = "WRF B",
                           "CC" = "WRF C"))

Seperate the data by gene target to ease layering in the final plot

#make three data layers
only_positives <<- subset(final_data, (!is.na(final_data$Facility)))
only_n1 <- subset(only_positives, target == "N1")
only_n2 <- subset(only_positives, target == "N2")
only_background <<-final_data %>% 
  select(c(date, cases_cum_clarke, new_cases_clarke, X7_day_ave_clarke)) %>%
  group_by(date) %>% summarise_if(is.numeric, mean)

#specify fun colors
background_color <- "#7570B3"
seven_day_ave_color <- "#E6AB02"
marker_colors <- c("N1" = '#1B9E77',"N2" ='#D95F02')

#remove facilty C for now
#only_n1 <- only_n1[!(only_n1$Facility == "WRF C"),]
#only_n2 <- only_n2[!(only_n2$Facility == "WRF C"),]

only_n1 <- only_n1[!(only_n1$Facility == "WRF A" & only_n1$date == "2020-11-02"), ]
only_n2 <- only_n2[!(only_n2$Facility == "WRF A" & only_n2$date == "2020-11-02"), ]

Build the main plot

      #first layer is the background epidemic curve
        p1 <- only_background %>%
              plotly::plot_ly() %>%
              plotly::add_trace(x = ~date, y = ~new_cases_clarke, 
                                type = "bar", 
                                hoverinfo = "text",
                                text = ~paste('</br> Date: ', date,
                                                     '</br> Daily Cases: ', new_cases_clarke),
                                alpha = 0.5,
                                name = "Daily Reported Cases",
                                color = background_color,
                                colors = background_color,
                                showlegend = FALSE) %>%
            layout(yaxis = list(title = "Clarke County Daily Cases", showline=TRUE)) %>%
            layout(legend = list(orientation = "h", x = 0.2, y = -0.3))
        
        #renders the main plot layer two as seven day moving average
        p1 <- p1 %>% plotly::add_trace(x = ~date, y = ~X7_day_ave_clarke, 
                             type = "scatter",
                             mode = "lines",
                             hoverinfo = "text",
                            text = ~paste('</br> Date: ', date,
                                                     '</br> Seven-Day Moving Average: ', X7_day_ave_clarke),
                             name = "Seven Day Moving Average Athens",
                             line = list(color = seven_day_ave_color),
                             showlegend = FALSE)
      

        
        #renders the main plot layer three as positive target hits
        
        p2 <- plotly::plot_ly() %>%
          plotly::add_trace(x = ~date, y = ~mean_copy_num_L,
                                       type = "scatter",
                                       mode = "markers",
                                       hoverinfo = "text",
                                       text = ~paste('</br> Date: ', date,
                                                     '</br> Facility: ', Facility,
                                                     '</br> Target: ', target,
                                                     '</br> Copies/L: ', round(mean_copy_num_L, digits = 2)),
                                       data = only_n1,
                                       symbol = ~Facility,
                                       marker = list(color = '#1B9E77', size = 8, opacity = 0.65),
                                       showlegend = FALSE) %>%
          plotly::add_trace(x = ~date, y = ~mean_copy_num_L,
                                       type = "scatter",
                                       mode = "markers",
                                       hoverinfo = "text",
                                       text = ~paste('</br> Date: ', date,
                                                     '</br> Facility: ', Facility,
                                                     '</br> Target: ', target,
                                                     '</br> Copies/L: ', round(mean_copy_num_L, digits = 2)),
                                       data = only_n2,
                                       symbol = ~Facility,
                                       marker = list(color = '#D95F02', size = 8, opacity = 0.65),
                                       showlegend = FALSE) %>%
            layout(yaxis = list(title = "SARS CoV-2 Copies/L", 
                                 showline = TRUE,
                                 type = "log",
                                 dtick = 1,
                                 automargin = TRUE)) %>%
            layout(legend = list(orientation = "h", x = 0.2, y = -0.3))
        
        #adds the limit of detection dashed line
        p2 <- p2 %>% plotly::add_segments(x = as.Date("2021-06-30"), 
                                          xend = ~max(date + 10), 
                                          y = 3571.429, yend = 3571.429,
                                          opacity = 0.35,
                                          line = list(color = "black", dash = "dash")) %>%
          layout(annotations = list(x = as.Date("2021-06-30"), y = 3.8, xref = "x", yref = "y", 
                                    text = "Limit of Detection", showarrow = FALSE))

        

        p1

        p2

Combine the two main plot pieces as a subplot

#seperate n1 and n2 frames by site
#n1
wrf_a_only_n1 <- subset(only_n1, Facility == "WRF A")
wrf_b_only_n1 <- subset(only_n1, Facility == "WRF B")
wrf_c_only_n1 <- subset(only_n1, Facility == "WRF C")

#n2
wrf_a_only_n2 <- subset(only_n2, Facility == "WRF A")
wrf_b_only_n2 <- subset(only_n2, Facility == "WRF B")
wrf_c_only_n2 <- subset(only_n2, Facility == "WRF C")


#rejoin the old data frames then seperate in to averages for each plant. 
wrfa_both <- full_join(wrf_a_only_n1, wrf_a_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))

## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke",
## "X7_day_ave_clarke", "Facility", "collection_num", "target",
## "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "se_L", "mean_total_copies",
## "sd_total_copies", "lo_95", "up_95", "log_copy_per_L")

wrfb_both <- full_join(wrf_b_only_n1, wrf_b_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))

## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke",
## "X7_day_ave_clarke", "Facility", "collection_num", "target",
## "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "se_L", "mean_total_copies",
## "sd_total_copies", "lo_95", "up_95", "log_copy_per_L")

wrfc_both <- full_join(wrf_c_only_n1, wrf_c_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))

## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke",
## "X7_day_ave_clarke", "Facility", "collection_num", "target",
## "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "se_L", "mean_total_copies",
## "sd_total_copies", "lo_95", "up_95", "log_copy_per_L")

#get max date
maxdate <- max(wrfa_both$date)
mindate <- min(wrfa_both$date)

Build loess smoothing figures figures

This makes the individual plots

#**************************************WRF A PLOT**********************************************
#add trendlines 
#extract data from geom_smooth
#both extract
# *********************************span 0.6***********************************
#*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_botha <- ggplot(wrfa_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_botha<<-..y..), method = "loess", color = '#1B9E77', 
              span = 0.25, n = 463)

## Warning: Ignoring unknown aesthetics: outfit

#look at the fits to align dates and total observations
#both
extract_botha

## `geom_smooth()` using formula 'y ~ x'

fit_botha

##   [1] 11.38411 11.43878 11.49245 11.54516 11.59693 11.64778 11.69774 11.74683
##   [9] 11.79508 11.84251 11.88915 11.93503 11.98016 12.02457 12.06829 12.11126
##  [17] 12.15340 12.19470 12.23513 12.27469 12.31335 12.35111 12.38794 12.42383
##  [25] 12.45877 12.49273 12.52571 12.55769 12.58864 12.61881 12.64838 12.67728
##  [33] 12.70546 12.73286 12.75940 12.78504 12.80970 12.83332 12.85584 12.87721
##  [41] 12.89734 12.91620 12.93370 12.95054 12.96732 12.98388 13.00005 13.01566
##  [49] 13.03054 13.04452 13.05742 13.06909 13.07935 13.08803 13.09496 13.09998
##  [57] 13.10291 13.10322 13.10068 13.09562 13.08832 13.07910 13.06826 13.05610
##  [65] 13.04293 13.02904 13.01476 13.00037 12.98619 12.97252 12.95967 12.94481
##  [73] 12.92539 12.90211 12.87564 12.84670 12.81596 12.78412 12.75188 12.71991
##  [81] 12.68892 12.65959 12.63262 12.60870 12.58851 12.56904 12.54708 12.52311
##  [89] 12.49760 12.47102 12.44383 12.41650 12.38950 12.36330 12.33836 12.31516
##  [97] 12.29416 12.27583 12.26063 12.24752 12.23513 12.22344 12.21245 12.20216
## [105] 12.19254 12.18360 12.17532 12.16770 12.16072 12.15438 12.14868 12.14359
## [113] 12.13911 12.13640 12.13636 12.13860 12.14273 12.14839 12.15518 12.16272
## [121] 12.17063 12.17853 12.18603 12.19274 12.19829 12.20230 12.20437 12.20512
## [129] 12.20547 12.20554 12.20543 12.20528 12.20518 12.20526 12.20563 12.20639
## [137] 12.20767 12.20958 12.21223 12.21574 12.22021 12.22532 12.23068 12.23630
## [145] 12.24222 12.24846 12.25504 12.26200 12.26935 12.27713 12.28536 12.29407
## [153] 12.30327 12.31301 12.32329 12.33416 12.34563 12.35773 12.37049 12.38393
## [161] 12.39999 12.42011 12.44361 12.46977 12.49789 12.52728 12.55724 12.58705
## [169] 12.61602 12.64345 12.66864 12.69089 12.70949 12.72374 12.73687 12.75235
## [177] 12.76983 12.78897 12.80943 12.83088 12.85296 12.87535 12.89769 12.91965
## [185] 12.94089 12.96106 12.97983 12.99685 13.01179 13.02430 13.03404 13.04068
## [193] 13.04386 13.04326 13.03852 13.02965 13.01712 13.00133 12.98265 12.96149
## [201] 12.93822 12.91324 12.88694 12.85970 12.83192 12.80397 12.77626 12.74916
## [209] 12.72308 12.69350 12.65646 12.61311 12.56461 12.51212 12.45680 12.39980
## [217] 12.34230 12.28545 12.23040 12.17833 12.13038 12.08771 12.05150 12.01745
## [225] 11.98090 11.94240 11.90253 11.86183 11.82087 11.78022 11.74043 11.70207
## [233] 11.66570 11.63188 11.60118 11.57415 11.55136 11.53093 11.51073 11.49093
## [241] 11.47173 11.45329 11.43581 11.41945 11.40439 11.39082 11.37892 11.36886
## [249] 11.36082 11.35499 11.35154 11.35182 11.35667 11.36547 11.37761 11.39247
## [257] 11.40946 11.42794 11.44732 11.46697 11.48630 11.50468 11.52151 11.53617
## [265] 11.54805 11.56017 11.57560 11.59386 11.61445 11.63687 11.66063 11.68523
## [273] 11.71018 11.73497 11.75912 11.78212 11.80349 11.82273 11.83933 11.85495
## [281] 11.87147 11.88875 11.90665 11.92503 11.94374 11.96265 11.98161 12.00048
## [289] 12.01913 12.03741 12.05517 12.07229 12.08861 12.10411 12.11889 12.13306
## [297] 12.14670 12.15990 12.17275 12.18534 12.19776 12.21011 12.22246 12.23491
## [305] 12.24756 12.26048 12.27377 12.28732 12.30092 12.31456 12.32822 12.34186
## [313] 12.35548 12.36904 12.38252 12.39591 12.40917 12.42229 12.43525 12.44801
## [321] 12.46057 12.47289 12.48495 12.49759 12.51150 12.52643 12.54215 12.55844
## [329] 12.57504 12.59172 12.60826 12.62440 12.63992 12.65458 12.66815 12.68037
## [337] 12.69103 12.70030 12.70859 12.71599 12.72261 12.72855 12.73391 12.73880
## [345] 12.74332 12.74757 12.75166 12.75567 12.75973 12.76393 12.76836 12.77280
## [353] 12.77694 12.78078 12.78433 12.78758 12.79055 12.79324 12.79565 12.79778
## [361] 12.79964 12.80123 12.80256 12.80363 12.80444 12.80454 12.80359 12.80173
## [369] 12.79912 12.79591 12.79224 12.78827 12.78416 12.78005 12.77609 12.77244
## [377] 12.76924 12.76666 12.76483 12.76243 12.75823 12.75259 12.74585 12.73835
## [385] 12.73043 12.72243 12.71470 12.70759 12.70142 12.69654 12.69331 12.69205
## [393] 12.69312 12.69654 12.70190 12.70888 12.71716 12.72640 12.73628 12.74649
## [401] 12.75668 12.76654 12.77574 12.78395 12.79086 12.79614 12.79945 12.80214
## [409] 12.80563 12.80972 12.81424 12.81900 12.82379 12.82844 12.83276 12.83655
## [417] 12.83963 12.84181 12.84290 12.84270 12.84104 12.83852 12.83586 12.83304
## [425] 12.83003 12.82682 12.82337 12.81967 12.81569 12.81141 12.80681 12.80186
## [433] 12.79653 12.79082 12.78468 12.77810 12.77106 12.76354 12.75550 12.74692
## [441] 12.73779 12.72822 12.71833 12.70809 12.69750 12.68652 12.67514 12.66334
## [449] 12.65110 12.63841 12.62524 12.61157 12.59738 12.58266 12.56731 12.55128
## [457] 12.53464 12.51741 12.49965 12.48140 12.46271 12.44363 12.42420

#assign fits to a vector
both_trenda <- fit_botha

#extract y min and max for each
limits_botha <- ggplot_build(extract_botha)$data

## `geom_smooth()` using formula 'y ~ x'

limits_botha <- as.data.frame(limits_botha)
both_ymina <- limits_botha$ymin
both_ymaxa <- limits_botha$ymax

#reassign dataframes (just to be safe)
work_botha <- wrfa_both

#fill in missing dates to smooth fits
work_botha <- work_botha %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_botha <- work_botha$date

#create a new smooth dataframe to layer
smooth_frame_botha <- data.frame(date_vec_botha, both_trenda, both_ymina, both_ymaxa)

#WRF A
#plot smooth frames
p_wrf_a <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_botha, y = ~both_trenda,
                    data = smooth_frame_botha,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_botha,
                                  '</br> Median Log Copies: ', round(both_trenda, digits = 2)),
                    line = list(color = '#1B9E77', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_botha, ymin = ~both_ymina, ymax = ~both_ymaxa,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_botha, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxa, digits = 2),
                                  '</br> Min Log Copies: ', round(both_ymina, digits = 2)),
                    name = "",
                    fillcolor = '#1B9E77',
                    line = list(color = '#1B9E77')) %>%
                layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF A") %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfa_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#1B9E77', size = 6, opacity = 0.65))

p_wrf_a

save(p_wrf_a, file = "./site_objects/wrf_a_year2.rda")

#**************************************WRF B PLOT**********************************************
#add trendlines 
#extract data from geom_smooth
#both extract
# *********************************span 0.6***********************************
#*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_bothb <- ggplot(wrfb_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_bothb<<-..y..), method = "loess", color = '#D95F02', 
              span = 0.25, n = 463)

## Warning: Ignoring unknown aesthetics: outfit

#look at the fits to align dates and total observations
#both
extract_bothb

## `geom_smooth()` using formula 'y ~ x'

fit_bothb

##   [1] 10.85583 10.93331 11.00951 11.08441 11.15797 11.23018 11.30099 11.37039
##   [9] 11.43835 11.50483 11.56982 11.63327 11.69517 11.75549 11.81419 11.87134
##  [17] 11.92702 11.98125 12.03406 12.08545 12.13545 12.18408 12.23135 12.27729
##  [25] 12.32192 12.36525 12.40730 12.44809 12.48765 12.52570 12.56203 12.59670
##  [33] 12.62979 12.66136 12.69149 12.72026 12.74773 12.77398 12.79907 12.82308
##  [41] 12.84609 12.86816 12.88936 12.90871 12.92535 12.93954 12.95152 12.96157
##  [49] 12.96992 12.97685 12.98260 12.98743 12.99161 12.99537 12.99899 13.00272
##  [57] 13.00681 13.01011 13.01141 13.01092 13.00883 13.00536 13.00069 12.99503
##  [65] 12.98858 12.98154 12.97411 12.96649 12.95889 12.95150 12.94452 12.93750
##  [73] 12.92986 12.92160 12.91274 12.90327 12.89321 12.88257 12.87135 12.85956
##  [81] 12.84721 12.83430 12.82085 12.80686 12.79234 12.77628 12.75790 12.73752
##  [89] 12.71546 12.69206 12.66764 12.64251 12.61702 12.59148 12.56622 12.54157
##  [97] 12.51784 12.49537 12.47449 12.45192 12.42475 12.39382 12.35999 12.32410
## [105] 12.28700 12.24954 12.21257 12.17693 12.14348 12.11306 12.08651 12.06469
## [113] 12.04845 12.03499 12.02108 12.00696 11.99286 11.97901 11.96562 11.95294
## [121] 11.94118 11.93057 11.92135 11.91373 11.90795 11.90423 11.90281 11.90369
## [129] 11.90665 11.91149 11.91805 11.92613 11.93556 11.94617 11.95777 11.97019
## [137] 11.98325 11.99676 12.01055 12.02445 12.03826 12.05401 12.07353 12.09645
## [145] 12.12237 12.15090 12.18165 12.21422 12.24822 12.28326 12.31894 12.35488
## [153] 12.39068 12.42595 12.46029 12.49332 12.52464 12.55386 12.58059 12.60443
## [161] 12.62813 12.65436 12.68262 12.71241 12.74325 12.77462 12.80604 12.83701
## [169] 12.86703 12.89560 12.92223 12.94643 12.96769 12.98552 13.00212 13.01985
## [177] 13.03850 13.05784 13.07765 13.09770 13.11778 13.13765 13.15710 13.17589
## [185] 13.19381 13.21063 13.22613 13.24009 13.25227 13.26246 13.27044 13.27597
## [193] 13.27884 13.27882 13.27569 13.26953 13.26072 13.24948 13.23603 13.22059
## [201] 13.20338 13.18463 13.16456 13.14338 13.12132 13.09860 13.07544 13.05207
## [209] 13.02869 13.00232 12.97034 12.93357 12.89283 12.84893 12.80268 12.75490
## [217] 12.70641 12.65801 12.61054 12.56479 12.52160 12.48176 12.44610 12.41077
## [225] 12.37179 12.32983 12.28554 12.23960 12.19265 12.14535 12.09837 12.05237
## [233] 12.00799 11.96591 11.92677 11.89125 11.85999 11.82961 11.79668 11.76180
## [241] 11.72557 11.68859 11.65148 11.61482 11.57923 11.54531 11.51365 11.48486
## [249] 11.45955 11.43832 11.42176 11.40934 11.39990 11.39316 11.38885 11.38667
## [257] 11.38635 11.38759 11.39012 11.39365 11.39789 11.40258 11.40741 11.41211
## [265] 11.41639 11.42220 11.43139 11.44354 11.45822 11.47501 11.49348 11.51320
## [273] 11.53376 11.55473 11.57567 11.59618 11.61581 11.63416 11.65078 11.66790
## [281] 11.68775 11.70995 11.73409 11.75979 11.78665 11.81429 11.84231 11.87032
## [289] 11.89793 11.92475 11.95038 11.97444 11.99653 12.01868 12.04293 12.06890
## [297] 12.09624 12.12456 12.15351 12.18271 12.21180 12.24041 12.26816 12.29470
## [305] 12.31965 12.34265 12.36332 12.38248 12.40122 12.41953 12.43745 12.45499
## [313] 12.47216 12.48897 12.50545 12.52161 12.53746 12.55302 12.56831 12.58334
## [321] 12.59812 12.61268 12.62703 12.64106 12.65468 12.66792 12.68079 12.69334
## [329] 12.70559 12.71757 12.72929 12.74080 12.75211 12.76326 12.77427 12.78517
## [337] 12.79599 12.80618 12.81528 12.82342 12.83075 12.83740 12.84352 12.84925
## [345] 12.85473 12.86009 12.86548 12.87104 12.87691 12.88322 12.89013 12.89741
## [353] 12.90474 12.91209 12.91944 12.92676 12.93401 12.94117 12.94822 12.95512
## [361] 12.96185 12.96837 12.97467 12.98071 12.98646 12.99186 12.99687 13.00154
## [369] 13.00593 13.01007 13.01403 13.01784 13.02156 13.02524 13.02892 13.03264
## [377] 13.03647 13.04044 13.04461 13.04891 13.05324 13.05757 13.06190 13.06620
## [385] 13.07045 13.07464 13.07874 13.08275 13.08665 13.09041 13.09402 13.09746
## [393] 13.10071 13.10443 13.10915 13.11465 13.12071 13.12712 13.13365 13.14010
## [401] 13.14624 13.15186 13.15674 13.16067 13.16342 13.16478 13.16453 13.16323
## [409] 13.16158 13.15956 13.15714 13.15433 13.15110 13.14745 13.14335 13.13879
## [417] 13.13375 13.12824 13.12221 13.11568 13.10861 13.10105 13.09303 13.08456
## [425] 13.07567 13.06635 13.05662 13.04648 13.03596 13.02505 13.01377 13.00213
## [433] 12.99014 12.97781 12.96515 12.95217 12.93889 12.92530 12.91142 12.89727
## [441] 12.88284 12.86812 12.85304 12.83762 12.82185 12.80572 12.78923 12.77238
## [449] 12.75516 12.73757 12.71960 12.70126 12.68253 12.66342 12.64395 12.62416
## [457] 12.60402 12.58355 12.56274 12.54157 12.52006 12.49818 12.47595

#assign fits to a vector
both_trendb <- fit_bothb

#extract y min and max for each
limits_bothb <- ggplot_build(extract_bothb)$data

## `geom_smooth()` using formula 'y ~ x'

limits_bothb <- as.data.frame(limits_bothb)
both_yminb <- limits_bothb$ymin
both_ymaxb <- limits_bothb$ymax

#reassign dataframes (just to be safe)
work_bothb <- wrfb_both

#fill in missing dates to smooth fits
work_bothb <- work_bothb %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_bothb <- work_bothb$date

#create a new smooth dataframe to layer
smooth_frame_bothb <- data.frame(date_vec_bothb, both_trendb, both_yminb, both_ymaxb)

#WRF B
#plot smooth frames
p_wrf_b <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_bothb, y = ~both_trendb,
                    data = smooth_frame_bothb,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothb,
                                  '</br> Median Log Copies: ', round(both_trendb, digits = 2)),
                    line = list(color = '#D95F02', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_bothb, ymin = ~both_yminb, ymax = ~both_ymaxb,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothb, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxb, digits = 2),
                                  '</br> Min Log Copies: ', round(both_yminb, digits = 2)),
                    name = "",
                    fillcolor = '#D95F02',
                    line = list(color = '#D95F02')) %>%
                layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF B") %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfb_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#D95F02', size = 6, opacity = 0.65))

p_wrf_b

save(p_wrf_b, file = "./site_objects/wrf_b_year2.rda")

#**************************************WRF C PLOT********************************************** #add trendlines #extract data from geom_smooth # *********************************span 0.6*********************************** #*****************Must always update the n = TOTAL NUMBER OF DAYS*************************

extract_bothc <- ggplot(wrfc_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_bothc<<-..y..), method = "loess", color = '#E7298A', 
              span = 0.25, n = 463)

## Warning: Ignoring unknown aesthetics: outfit

#look at the fits to align dates and total observations
#both
extract_bothc

## `geom_smooth()` using formula 'y ~ x'

fit_bothc

##   [1] 10.86751 10.92648 10.98447 11.04148 11.09749 11.15251 11.20654 11.25956
##   [9] 11.31157 11.36258 11.41257 11.46154 11.50948 11.55640 11.60228 11.64716
##  [17] 11.69108 11.73402 11.77598 11.81695 11.85694 11.89593 11.93392 11.97091
##  [25] 12.00689 12.04186 12.07581 12.10874 12.14064 12.17142 12.20103 12.22950
##  [33] 12.25685 12.28311 12.30832 12.33251 12.35570 12.37794 12.39924 12.41964
##  [41] 12.43917 12.45787 12.47575 12.49268 12.50850 12.52325 12.53696 12.54967
##  [49] 12.56142 12.57223 12.58214 12.59118 12.59939 12.60681 12.61346 12.61938
##  [57] 12.62461 12.62885 12.63181 12.63358 12.63422 12.63382 12.63245 12.63018
##  [65] 12.62708 12.62323 12.61870 12.61357 12.60792 12.60180 12.59531 12.58675
##  [73] 12.57473 12.55976 12.54235 12.52301 12.50224 12.48055 12.45845 12.43646
##  [81] 12.41508 12.39481 12.37618 12.35968 12.34583 12.33254 12.31759 12.30126
##  [89] 12.28384 12.26561 12.24687 12.22788 12.20895 12.19036 12.17238 12.15531
##  [97] 12.13944 12.12504 12.11241 12.10071 12.08896 12.07723 12.06556 12.05401
## [105] 12.04264 12.03152 12.02068 12.01020 12.00013 11.99052 11.98144 11.97294
## [113] 11.96507 11.95811 11.95220 11.94721 11.94298 11.93938 11.93629 11.93354
## [121] 11.93102 11.92858 11.92607 11.92338 11.92034 11.91683 11.91271 11.90758
## [129] 11.90131 11.89415 11.88636 11.87820 11.86993 11.86182 11.85411 11.84707
## [137] 11.84096 11.83603 11.83255 11.83078 11.83097 11.83205 11.83284 11.83346
## [145] 11.83403 11.83466 11.83546 11.83655 11.83804 11.84004 11.84267 11.84605
## [153] 11.85029 11.85550 11.86179 11.86929 11.87810 11.88834 11.90012 11.91356
## [161] 11.93047 11.95211 11.97774 12.00662 12.03804 12.07125 12.10553 12.14014
## [169] 12.17435 12.20743 12.23865 12.26727 12.29257 12.31382 12.33483 12.35963
## [177] 12.38774 12.41866 12.45193 12.48707 12.52360 12.56104 12.59891 12.63674
## [185] 12.67404 12.71035 12.74517 12.77804 12.80846 12.83598 12.86010 12.88036
## [193] 12.89626 12.90734 12.91311 12.91439 12.91248 12.90770 12.90032 12.89064
## [201] 12.87896 12.86557 12.85076 12.83482 12.81804 12.80073 12.78317 12.76565
## [209] 12.74847 12.72726 12.69827 12.66270 12.62173 12.57658 12.52845 12.47852
## [217] 12.42800 12.37810 12.33000 12.28491 12.24403 12.20856 12.17969 12.15344
## [225] 12.12528 12.09560 12.06483 12.03335 12.00156 11.96988 11.93870 11.90842
## [233] 11.87944 11.85218 11.82702 11.80438 11.78465 11.76704 11.75045 11.73484
## [241] 11.72016 11.70639 11.69346 11.68135 11.67001 11.65940 11.64947 11.64018
## [249] 11.63149 11.62336 11.61575 11.60924 11.60437 11.60099 11.59897 11.59817
## [257] 11.59845 11.59968 11.60171 11.60441 11.60765 11.61128 11.61516 11.61917
## [265] 11.62315 11.62820 11.63532 11.64422 11.65463 11.66627 11.67885 11.69211
## [273] 11.70575 11.71951 11.73311 11.74626 11.75869 11.77012 11.78027 11.78987
## [281] 11.79982 11.81007 11.82060 11.83136 11.84232 11.85345 11.86471 11.87606
## [289] 11.88746 11.89889 11.91030 11.92166 11.93294 11.94517 11.95925 11.97489
## [297] 11.99178 12.00965 12.02821 12.04715 12.06620 12.08507 12.10345 12.12107
## [305] 12.13764 12.15285 12.16643 12.17951 12.19330 12.20767 12.22246 12.23752
## [313] 12.25271 12.26787 12.28287 12.29754 12.31174 12.32533 12.33815 12.35005
## [321] 12.36089 12.37052 12.37878 12.38572 12.39155 12.39639 12.40039 12.40365
## [329] 12.40633 12.40854 12.41041 12.41207 12.41365 12.41528 12.41709 12.41921
## [337] 12.42176 12.42431 12.42638 12.42803 12.42932 12.43030 12.43104 12.43159
## [345] 12.43201 12.43237 12.43271 12.43311 12.43362 12.43430 12.43520 12.43546
## [353] 12.43432 12.43199 12.42869 12.42464 12.42007 12.41520 12.41023 12.40540
## [361] 12.40092 12.39701 12.39390 12.39179 12.39092 12.39043 12.38938 12.38792
## [369] 12.38616 12.38424 12.38228 12.38040 12.37874 12.37742 12.37657 12.37630
## [377] 12.37676 12.37807 12.38034 12.38386 12.38868 12.39460 12.40145 12.40902
## [385] 12.41713 12.42559 12.43422 12.44281 12.45118 12.45915 12.46651 12.47309
## [393] 12.47869 12.48462 12.49215 12.50097 12.51077 12.52127 12.53217 12.54317
## [401] 12.55397 12.56427 12.57378 12.58220 12.58922 12.59456 12.59792 12.60010
## [409] 12.60209 12.60387 12.60540 12.60666 12.60761 12.60823 12.60848 12.60834
## [417] 12.60778 12.60676 12.60525 12.60324 12.60068 12.59765 12.59424 12.59045
## [425] 12.58629 12.58175 12.57683 12.57154 12.56588 12.55985 12.55345 12.54668
## [433] 12.53954 12.53203 12.52416 12.51592 12.50732 12.49836 12.48904 12.47935
## [441] 12.46931 12.45890 12.44811 12.43694 12.42540 12.41348 12.40119 12.38852
## [449] 12.37549 12.36209 12.34832 12.33419 12.31969 12.30483 12.28962 12.27409
## [457] 12.25822 12.24201 12.22546 12.20855 12.19129 12.17365 12.15565

#assign fits to a vector
both_trendc <- fit_bothc

#extract y min and max for each
limits_bothc <- ggplot_build(extract_bothc)$data

## `geom_smooth()` using formula 'y ~ x'

limits_bothc <- as.data.frame(limits_bothc)
both_yminc <- limits_bothc$ymin
both_ymaxc <- limits_bothc$ymax

#reassign dataframes (just to be safe)
work_bothc <- wrfc_both

#fill in missing dates to smooth fits
work_bothc <- work_bothc %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_bothc <- work_bothc$date

#create a new smooth dataframe to layer
smooth_frame_bothc <- data.frame(date_vec_bothc, both_trendc, both_yminc, both_ymaxc)

#WRF C
#plot smooth frames
p_wrf_c <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_bothc, y = ~both_trendc,
                    data = smooth_frame_bothc,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothc,
                                  '</br> Median Log Copies: ', round(both_trendc, digits = 2)),
                    line = list(color = '#E7298A', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_bothc, ymin = ~both_yminc, ymax = ~both_ymaxc,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothc, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxc, digits = 2),
                                  '</br> Min Log Copies: ', round(both_yminc, digits = 2)),
                    name = "",
                    fillcolor = '#E7298A',
                    line = list(color = '#E7298A')) %>%
                layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF C") %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfc_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#E7298A', size = 6, opacity = 0.65))

p_wrf_c

save(p_wrf_c, file = "./site_objects/wrf_c_year2.rda")

keeping in case

#save(wrfa_both, file = "./plotly_objs/wrfa_both.rda")
#save(wrfb_both, file = "./plotly_objs/wrfb_both.rda")
#save(wrfc_both, file = "./plotly_objs/wrfc_both.rda")

#save(date_vec_botha, file = "./plotly_objs/date_vec_botha.rda")
#save(date_vec_bothb, file = "./plotly_objs/date_vec_bothb.rda")
#save(date_vec_bothc, file = "./plotly_objs/date_vec_bothc.rda")

#save(both_ymina, file = "./plotly_objs/both_ymina.rda")
#save(both_ymaxa, file = "./plotly_objs/both_ymaxa.rda")

#save(both_yminb, file = "./plotly_objs/both_yminb.rda")
#save(both_ymaxb, file = "./plotly_objs/both_ymaxb.rda")

#save(both_yminc, file = "./plotly_objs/both_yminc.rda")
#save(both_ymaxc, file = "./plotly_objs/both_ymaxc.rda")